Transitions in non-conserving models of Self- Organized Criticality 



Stefano Lise a and Henrik Jeldtoft Jensen b 
Department of Mathematics, Imperial College 
180 Queen's Gate, London SW7 2BZ 
United Kingdom 

Abstract 

We investigate a random-neighbours version of the two dimensional non- 
conserving earthquake model of Olami, Feder and Christensen [Phys. Rev. 
Lett. 68, 1244 (1992)]. We show both analytically and numerically that 
criticality can be expected even in the presence of dissipation. As the critical 
level of conservation, a c , is approached, the cut-off of the avalanche size 
distribution scales as £ ~ (a c — a) -3 / 2 . The transition from non-SOC to SOC 
behaviour is controlled by the average branching ratio a of an avalanche, which 
can thus be regarded as an order parameter of the system. The relevance of 
the results are discussed in connection to the nearest-neighbours OFC model 
(in particular we analyse the relevance of synchronization in the latter). 
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In what is now a commonly referenced paper |IJ, Bak, Tang and Wiesenfeld proposed that 
extensive systems of many coupled elements naturally evolve toward a dynamical critical 
state. They named this phenomenon Self-Organized Criticality and tested the idea on a 
cellular automata, the so called sandpile model. It was shown in that paper, and many 
others which followed J2], that, independently of the starting, configuration the sandpile 
model evolves to a stationary state which does not possess any characteristic spatial scale. 
The amplitude of the response of the system to an external perturbation follows a power 
law distribution. 

From a theoretical point of view the main effort has been directed towards understanding 
which are the fundamental mechanisms that lead a system to be self-organized critical. It is 
well known, for instance, that for a particular class of sandpile models a necessary condition 
is a conserving local dynamics ||(f|: any amount of dissipation introduces a characteristic 
length-scale and thus destroys criticality. Nonetheless this does not seem to be a universal 
requirement, since in some other classes of models, conservation does not appear to be 
needed. In this respect the earthquake model introduced by Olami, Feder and Christensen || 
is of particular interest. In this model it is possible to directly control the level of conservation 
of the dynamics, through a parameter a. When a = 1/4 the system is conserving and it 
probably belongs to the same universality class as the sandpile model jy. But, contrary 
to the latter, there are converging evidences that the OFC model remains critical even 
when dissipation is introduced (a < 1/4) |5HlO|. By reducing the value of a one should 



be able in principle to observe a transition from S.O.C. behaviour to non-S.O.C, although 
a general agreement on the critical value , a c , does not yet exist. The spectrum of values 



ranges from the recently proposed a c = [10] to a c ~ 0.18 Hi, passing through the 



originally estimated a c ~ 0.05 [|]|| fl2| . It has been suggested that the observed criticality 
in this model is a consequence of the "imperfect" synchronization in the system |ID| , |II| . The 
strong correlations so induced would allow an avalanche (i.e. the response to an external 
perturbation) to be infinitely big, although it keeps on losing "energy" . 

In this letter we investigate the random neighbours version of the OFC model ||. We 



find that a sharp transition from non-S.O.C. to S.O.C. behaviour occurs at a c ~ 2/9 in this 
model. Even more interesting, this transition is correctly described by the branching ratio, 
which can thus be regarded as an order parameter of the system, in analogy to equilibrium 
phase transitions. The results obtained for the random neighbours model are interesting in 
their own right, but their relevance is fully recognised when they are related to the nearest 
neighbours version. In fact they provide not only a mean-field description of the OFC model, 
but also properly address the role played by synchronization. For a system where the sole 
mechanism leading to SOC is synchronization, one would expect the random neighbours 
version to be non critical, since the latter cannot synchronize. 

The OFC model is a a coupled-map lattice model which, despite its simplicity, is thought 
to capture some of the essential features of earthquake dynamics ||13|| . To each site of a 2- 



dimensional square lattice is associated a real continuous "energy" E{. The system is driven 
continuously and uniformly, which means that all the Ei values are simultaneously increased 
with time at the same rate. Avalanche (or earthquake) dynamics is simulated by assuming 
that a single site is unable to store more than a finite amount of energy E c . As soon as a 
site becomes unstable (i.e. Ei > E c ) an avalanche is triggered: the global driving is stopped 
and the system evolves according to the following local relaxation rule: 

I Ei -> 

if Ei > E c => I (1) 
I E nn — > E nn + aiEi 

until all of sites are below E c . In eq.fll]) "nn" stands for the collection of nearest neighbours 
to site i. The parameter a G [0, 4] controls the conservation level of the dynamics (a = ~ 
corresponds to the conservative case). The random neighbours version we consider differs 
from the OFC model only in the choice of neighbours: an unstable site distributes an energy 
aEi to 4 randomly chosen sites. It is well known that boundary conditions play a crucial 
role in SOC |TJ]]. In accordance with previous studies we use open boundary conditions so 
that if one of the random neighbours is a boundary site the energy aEi is simply lost. 

The presence of a transition from localised to non-localised behaviour of the avalanches 
can be understood through the following argument. For notational clarity in the following 



we will distinguish between a stable and an unstable site through the superscripts — and 
+ . Let P + (E + ) be the probability that a random site will become active as a consequence 
of receiving a contribution of magnitude aE + . P + (E + ) is equal to the probability of a site 
having a value between E c — aE + and E c . The average number of active sites produced by 
un unstable site with energy E + is then 4P + (E + ). The branching ratio a is defined as the 
average number of new active sites created by an unstable site. Averaging over the whole 
spectrum of possible P + (E + ) (that is an average over E + ) we get for a 



where P(E) is the distribution of the dynamical variable in the system. Obviously the 
chance to have an infinite avalanche is directly related to the condition a > 1. 
An exact analytic calculation of P + is, unfortunately, too complicated, since it involves a 
detailed knowledge of P{E) (see fig.l). In order to continue we have to assume a specific 
functional form for P{E). For simplicity we approximate the distribution of subcritical E- 
values by a uniform distribution on the interval [0, E c ]. We obtain P + (E + ) = aE + /E c and 
therefore P + = a(E + )/E c where (E + ) is the average value of a collapsing site. Consequently 
it follows that 



It is immediately clear from eq. (|3]) that a will be greater then 1 even for some values of 
a < 1/4 since (E + ) > E c . The condition a > 1 can also be read as a conservation law. 
It states that the average contribution 4a(E + ) to the avalanche from the collapse of an 
unstable site must be greater or equal to E c otherwise the avalanche will die exponentially. 

We now estimate (E + ). Consider an active site i and a site j which will be active as a 
consequence of the action of i. We have Ef = Ej + aEf , where the superscripts + and — 
in the site j distinguish between two successive time steps. Taking the average on both side 
of this equation we get (Ej~) = {Ej) +a{Ef). By assumption E c — aEf < Ej < E c . Thus, 
if we again approximate the distribution of ^-values of the subcritical sites by a uniform 
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(2) 



a = 4a 



(£+> 



(3) 



E c 



distribution on [0, E c ] we will have (Ej ) = E c — ja(E^). Moreover, since E^ is a new active 
site (Ej~) = (Ef). Combining these expressions we obtain 

(E + ) = ^ (4) 

2 

The branching ratio is accordingly given by 

° = Y^ (5) 

1 2 

The condition for infinite avalanches a > 1 is then 

2 

a > a c = - ~ 0.222... (6) 
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The above calculation, although very crude in the assumption of a uniform distribution of 
subcritical .E-values, predicts two different phases for the random neighbours model: a "low- 
ex" phase, where avalanches are essentially smaller than a characteristic size, and a "high- 
er" phase, where on the contrary it is possible to have avalanches of any sizes. Computer 
simulations indeed confirm this point. Figure 2 shows the avalanche size distribution for 
increasing values of a. The size is measured by counting the total number of toppling in 
one avalanche. For a < 0.22 the distribution are very well described by the function 

P Q ( S )oc S ~ 3 / 2 exp(^) (7) 

where the cut-off £, for sufficently large system, is independent of L (L is linear dimension 
of the system). Fig. 3 shows the fit of the simulated P a (s) distributions to the exponential 
form in Eq. [7[ For a > 0.23 the cut-off in P a (s) scales with L. This is the signature 
of criticality. In fig. 4 we show £ as a function of a. The data fits well the expression 
£ ~ (a c — a) -1 ' 5 , with a c = 0.2255. Note that £ ~ (1/4 — a)~ 2 was found in the study of 



the random neighbour sand pile model fT5 |, i.e. the coherence length is only infinite in this 
model when conservation is established at a = 1/4. Even more exciting is the behaviour of 
the measured branching ratio. In fig. 5 we report the measured branching ratio for different 
a as a function of 1/L. One sees that for a < 0.22 the graphs extrapolate to a value lim a 
well below 1. As soon as a > 0.23 lim a ~ 1 in accordance with Eq. || 

L^oo 



In conclusion, by identifying a transition at a finite conservation level we have shown 
that a non- conserving system can indeed be critical. The agreement between the analytical 
calculation and computer simulation is remarkable and somewhat surprising. The assump- 
tion of a uniform distribution of subcritical .E-values is not a very accurate approximation, 
see Fig.l. In this connection it is interesting to note the observation made by Pietronero, 
Tartaglia, and Zhang |L|]. These authors found that the form of energy distribution de- 
pends very much on the energy partition rule (i.e. on how the energy of an unstable site is 
distributed to its neighbours), whereas the average energy in the system is quite universal. 
This is a very interesting point, which would deserve further investigations and which might 
also explain why the uniform approximation for P(E) gives such an accurate estimate of 
a c . It might be possible in fact that for a particular energy partition rule the approximation 
becomes almost correct (or even exact). In [|16| was shown, for instance, that a random 
partition makes the peaks of Fig. 1 disappear. We believe this property of universality 
might be related to the very nature of S.O.C: measurable quantities should be independent 
of the details of the model. We have shown that synchronization is not the only mechanism 
present. Rather, the criticality appears to be related to the dynamics. This of course, does 
not prevent synchronization to be relevant for lower values of a as has been suggested for 



the nearest neighbour version of the OFC model ||T0| , |TT|| . Finally we have shown that the 
branching ratio plays the role of an order parameter in the considered model. We believe 
that this last point might also be useful for the analysis of other models. The natural step 
is to use it for the nearest neighbour version the OFC model. Preliminary results are en- 
couraging even though simulations show that the scaling properties of the branching ratio 
are more complicated than is the case in the random neighbour model. 
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Captions 

Figure 1. 

Energy distribution per site, P(E), for a = 0.23 and L = 400 (continuous line) (E c = 1). 
The step function (dashed line) is the approximation used in the calculations. At E — 0, 
P{E) extends up to ~ 33. 
Figure 2. 

Avalanche size distribution for L = 100, 200 and 400 and for (a) a = 0.20, (b) a = 0.21, 
(c) a = 0.22 and (d) a = 0.23. 
Figure 3. 

Avalanche size distributions for a < a c . The ditributions do not scale with system size. The 
continuous lines are fits of the form P a (s) oc s~ 3 / 2 exp(|) 
Figure 4. 

Cutoff in avalanche size distribution, £, as a function of a c — a (a c = 0.2255). The solid line 
is the interpolation £ ~ (0.2255 - a)~ L5 
Figure 5. 

Measured branching ratio as a function of 1/L (L — 100, 200, 400), for different values of 
a. From bottom to top, a = 0.2,0.21,0.22,0.225,0.23,0.24,0.25. One can clearly see the 
abrupt change in behaviour for a ~ 0.225. 
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